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In the presence of fermionic matter the topologically distinct vacua of the 
^ . standard model are metastable and can decay by tunneling through the sphaleron 

■ barrier. This process annihilates one fermion per doublet due to the anomalous 

non-conservation of baryon and lepton currents and is accompanied by a pro- 
duction of gauge and Higgs bosons. We present a numerical method to obtain 
local bounce solutions which minimize the Euclidean action in the space of all 
configurations connecting two adjacent topological sectors. These solutions de- 
termine the decay rate and the configuration of the fields after the tunneling. 
We also follow the real time evolution of this configuration and analyze the spec- 
trum of the created bosons. If the matter density exceeds some critical value, 
the exponentially suppressed tunneling triggers off an avalanche producing an 
enormous amount of bosons. 



1 diakonov@lnpi.spb.su 

2 goeke@hadron.tp2.ruhr-uni-bochum.de 



1 Introduction 



Baryon and lepton number violating processes in the electroweak theory have 
been the subject of many recent investigations. They are due to the anomaly of 
the baryon and lepton currents, discovered by 't Hooft [U, and the non-trivial 
topological structure of the electroweak theory. Faddeev Q and Jackiw and 
Rebbi Q found that the potential energy is periodic in a certain functional of 
the field, the Chern-Simons number Nqs, so that instead of one unique vac- 
uum there exist infinitely many field configurations with zero potential energy, 
classified by integer values of Nqs ■ 

Each transition between vacua with AiVcs = 1 is accompanied by a change 
of the baryon and lepton number by one unit per fermion generation. The vacua 
are separated by an energy barrier, called sphaleron barrier [|], |H| , whose height 
is of the order of 10 TeV. Under ordinary conditions the barrier can only be 
overcome by tunneling, but the tunneling probability is suppressed by the factor 
exp(— 25i nst ) ~ 10~ 153 with the instanton action 5i nst = 87T 2 /g 2 , g rj 0.67, 
which means that the process practically never happens. 

Under special conditions, however, the fermion number violation rate might 
well be significant. For example, a large temperature (of the order of mw) 
allows the system to cross the barrier classically || 0, |j ; this process might have 
played a key role for the generation and conservation of the baryon asymmetry 
in the early universe. The energy which is necessary to overcome the barrier 
can possibly also be provided by the incoming particles in a collision if the 
particle energy is of the order of 10 ~ 100 TeV || [l(], 11, 12]. Hence, fermion 
number violation might be observable at future supercolliders. 

In this paper we will investigate a third possibility to obtain fermion number 
violation at a reasonable rate, namely in a surrounding of high density matter 



14]. The mechanism of how the suppression of the transition rate is reduced 



is as follows: Matter of high density is described by a chemical potential /i, 
which is, at temperature zero, the energy up to which the Fermi-levels are 
filled. A transition with ANqs = 1 creates a fermion which has to be placed 
into the first free level, i.e. it has energy [i. This energy must be added to the 



potential energy of the gauge and Higgs bosons [|TJ, 14]: 

<ot = ^pot+^V C S • (1-1) 

The extra term causes the previously degenerate vacua to become metastable 
and the height of the barrier to be reduced (see Fig. 1 where V^ ot is plotted for 
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a negative fx) so that the tunneling probability increases. At a certain critical 
value /x cr it, the barrier disappears completely and the states at integer Ncs 
become unstable. 




Fig. 1: Schematic plot of the tunneling process between two topological 
sectors. 

The decay rate of a metastable state per volume can be calculated by the 
semi-classical WKB method. Following Coleman it is expressed in the 
form T/V = B exp(-A). In order to find the exponent A one has to solve 
the classical Euclidean equations of motion, i.e. one has to find the classical 
motion of the system in the potential —Vp Qt from some metastable ground 
state (e.g. a state with Nqs = 0, as indicated in Fig. 1) to a configuration with 
the same potential energy which is on the other side of the barrier. We call 
this configuration "escape point". The exponent A is twice the action of this 
motion, minimized over all possible escape points. In Euclidean space, after 
reaching the escape point the system would move back the same way to the 
original ground state, so that this process is called "the bounce" |l5|]. Hence, 
the bounce represents the minimum of the Euclidean action in the space of all 
possible paths from a metastable ground state to the other side of the barrier. 
The prefactor B basically contains zero mode factors and the determinant of 
small fluctuations about the bounce, but in this work we will not be concerned 
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with it and only compute the exponent A. 

It turns out that the chemical potential (and hence the matter density) 
which is needed to obtain a reasonable decay rate is quite large, so that it can 
hardly occur under normal circumstances. Nevertheless, the results may still 
be relevant because the decay rate is closely related to the rate of baryon and 



lepton number violation at high particle energies [14]. Also, to our knowledge, 
this is the first investigation of the metastable vacuum decay in a real theory, 
beyond the so-called thin wall approximation; and since it is a hard technical 
problem, our numerical procedure might be useful in other applications, ranging 
from spontaneous decay of heavy nuclei to inflationary scenarios of the early 
universe. 

In principle there are two different possibilities to find the bounce numeri- 
cally: Either one solves the equations of motion by some initial value method 
like Runge-Kutta, or one considers the action as a functional of the fields and 
minimizes it in the space of all fields configurations. The first method seems to 
be unfeasible because the potential —Vp Qt has no lower bound so that a slight 
deviation from the correct path will cause the system to fall into some abyss 
of the potential — Vp Qt . We have therefore decided to take the second route. 
For the numerical implementation of the minimization we use a discretization 



in space and time, based on a procedure presented by Adler and Piran [16], so 
that the action becomes a function of the values of the boson fields at the grid 
nodes. At each step one considers the action as a function of a certain field at 
a given point and keeps the other values fixed. A single step of a Newtonian 
algorithm is performed, then one takes another field or moves to the next point 
until one has performed a "sweep" through the whole lattice and starts with 
the next one. 

Related problems have already been treated in the literature, for example 
in [|l7]] the decay rate of a metastable state was also calculated by minimizing 
the Euclidean action on a grid in the context of technibaryons in the Skyrme 
model. In |l8| the Euclidean action of the electroweak theory is minimized, but 
with respect to only a few parameters in the space of parametrized functions. A 
method how to minimize the action of a model with scalar fields was presented 



in ||19|| . A short letter about the present method has been published recently 
& 

Another aim of this work is to investigate the fate of the system after the 



tunneling process has happened. It is known [15| that the motion in Euclidean 
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space does not only provide the probability of the barrier penetration, but also 
yields the most probable field configuration in which we will find the system 
after the tunneling. This is just the escape point of the bounce trajectory, 
i.e. the configuration at the other side of the barrier which belongs to the path 
in Euclidean space with the least action. The potential energy of the system 
after the tunneling is greater than the energy of the metastable minimum in its 
current topological sector (see Fig. 1). The difference \fi\ is the energy of the 
annihilated fermion which is now at the disposal of the boson fields. 

The field configuration then performs a real time motion in Minkowski space. 
Knowing the escape point field configuration and the field velocities (which are 
zero), one has to solve the standard Cauchi problem, which we perform by 
straightforward integration of the second-order differential equations. Usually 
the system will fall towards the minimum in the sector of the escape point, 
and eventually settles at that minimum. The original energy |^| of fermions 
from the Fermi surface is converted into classical radiation of bosons. In 3+1 
dimensions the amplitude of the outgoing wave package falls off as 1 jt so that 
the generally nonlinear equations of motion can be linearized at large times t, 
and the outgoing fields take the form of spherical waves. In order to evaluate 
the particle content of the multi-boson final state one has to carry out a Fourier 
transformation of the fields after their amplitudes got small enough. At this 
point we basically follow the work of Hellmund and Kripfganz |]2l|| (see also 
]), who took a slightly disturbed sphaleron as starting configuration, let it 
evolve in real time and analyzed the resulting multi-boson state. Our main 
modifications are the replacement of the sphaleron by the escape point of the 
bounce and the introduction of the chemical potential according to eq. (|1.1| ). 

In certain cases, however, when the fermion density is large enough though 
less than the critical density fi crit , the real time evolution of the fields is utterly 
different. The fields will not quiet down at the minimum but "splash" over 
the next barrier leading to a still lower minimum (at Nqs = 2 in Fig. 1), 
and so forth. As a result the fermion sea will be completely "dried out", and 
a huge amount of energy will be released in form of boson radiation. Naively, 
one would think that such an avalanche happens only when the fermion density 
exceeds the critical // crit when the system is allowed to roll down classically. It is 
amusing that actually the avalanche-like fermion annihilation can be triggered 
off by a spontaneous tunneling process at the first stage. 

The paper is organized as follows: In Section 2 we set up the model and 
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describe how matter of high density allows baryon and lepton number violating 
processes. In Section 3 we present in detail our numerical procedure to find the 
bounce trajectory. The classical motion after the tunneling is investigated in 
Section 4. In Section 5 we give the numerical results of our calculations, and 
finally we summarize our work in Section 6. 

2 Decay of high— density matter by barrier 
penetration 

We consider the minimal version of the standard electroweak theory with one 
Higgs doublet in the limit of vanishing Weinberg angle. We work with dimen- 
sionless rescaled quantities, the corresponding physical quantities can in general 
be obtained by multiplication with appropriate powers of the gauge boson mass 
mw- Sometimes this factor is already included in the definition, for details see 
||. The bosonic part of the Lagrangian is 

C = r ^f(^-\F^F a ^ + §(A**)t - ±u 2 (&$ _ 4 ) 2 ) (2.1) 

with the covariant derivative = — iA^ , = ^A^T a , the field strength 
Fp> = W = i[D„,D v \ , F£ v = d^Al - d v Al + ^ C A\A%, and the Higgs 
doublet $ = (*o)- v = run /my/ is the ratio of Higgs and gauge boson masses; 
the T a are the Pauli matrices. We work entirely in temporal gauge, Aq = 0, 
which restricts possible gauge transformations to time-independent ones, 

Ai -> U{Ai + idi)U ] , $ -> U<S> with U = U{r) G SU{2) . (2.2) 

The potential and kinetic energy and the Chern-Simons index are 

"1(^)2 + l( A $)t (A$) + - Af 

~\(A<tf + \&$\ , (2.3) 

Ncs = it* I d " r \- €ijk ( A ' djAt + h^^l 



r m W f ,3 

Tkin = —^T jar 



Ncs is only well-defined if the configuration space can be identified with the 
sphere S3, which requires the fields to be continuous at infinity. We will always 
fix them to the trivial vacuum there (Ai = 0, $ = (2)). For vacuum configu- 
rations Ai = iUdiU^ with U G 577(2) (pure gauge), Ncs is the integer winding 
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number of the mapping S3 — > SU(2) ~ S3. Adjacent topological sectors are 
separated by energy barriers. The field configurations which minimize the en- 
ergy Vp t for given Nqs have been calculated numerically by Akiba, Kikuchi, 
and Yanagida (AKY) ||. 

The fermions are coupled to the gauge and Higgs fields via the covariant 
derivative and the Yukawa coupling, respectively. We do not consider them 
explicitly, but note that due to the anomaly of the fermionic currents their 
number is not conserved, it varies with the Chern-Simons number of the clas- 
sical bosonic fields. For each doublet (i) the fermion number changes as 

ANi = AN CS • (2.4) 

We assume to have a macroscopic amount of fermions of a very high density in 
thermal equilibrium at zero (or very low) temperature. There are N g (N c + 1) = 
12 doublets, three leptonic and nine quark ones. We describe them by the 
chemical potentials /ij of the doublets; since we have zero temperature, fii is 
the energy up to which the fermionic levels are occupied (Fermi energy). In a 
process connecting two adjacent vacua one level crosses the gap, all others are 
shifted such that after the transition we have the same spectrum again. But 
now one more level is occupied or depleted, depending on whether the levels 
went up or down, hence the created or annihilated fermion has the energy 
Eq. ( |2.4| ) allows to include the change of the fermionic energy into the bosonic 
energy functional. The new potential energy is 14] 



Cot = ^ P ot + /xiVcs (2.5) 



with Vp 0t and iV C g from (2^3). [i = J2i Hi 1S the sum of the chemical potentials 
of the single doublets. Obviously, we have fixed the zero-point of the energy to 
the trivial vacuum with Nqs = . Besides, we can neglect the change of the 
Fermi energy due to the creation or annihilation of a small number of fermions. 

The additional term modifies the curvature of the potential around the 
ground states of the topological sectors. A calculation of the second derivative 
of Vp Qt with respect to the fields yields the following modes in momentum space 



2 

1 + p 2 , v 2 + p 2 , and 1 + p 2 ± -4r^ |p| ; (2.6) 
besides, there is a zero mode due to gauge freedom. We see that there are 
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negative modes if \fi\ exceeds the critical value 

16vr 2 . „. 

Merit = — r m w • ( 2 - 7 ) 

Therefore, at > // crit the system may roll down classically without any 



tunneling. A numerical calculation [24] shows that indeed the minimal energy 
barriers between the topological sectors vanish in this case. 

We are mainly interested in the case < \fi\ < /v it , where the topologically 
distinct minima with integer Nqs have different energies Vp Qt = ^Nqs and 
are separated by energy barriers. Therefore, they are only local minima and 
are thus metastable and can decay spontaneously by quantum tunneling to 
the adjacent topological sector with lower ground state energy. The energy 
difference between the two ground states is the energy of the bosonic fields 
after the tunneling that will eventually be carried away by the outgoing boson 
radiation. 

As mentioned in the introduction, the tunneling rate per volume is of the 
form 1 15] 

T/V = Bexp(-A) (2.8) 

where A is twice the Euclidean action of the bounce trajectory, which is the 
classical path that connects the decaying state and the turning point at the 
other side of the inverted barrier and minimizes this action. The prefactor B 
can be found from the small oscillation determinant about the bounce, with 
one negative mode and the zero modes removed. B also contains the Jacobian 
factors of the transformation groups which leave the action invariant, i.e. which 
correspond to the zero modes. The factors coming from translational invariance 
in space and time are not included in B but in the left hand side of ( |2.8j ), they 
lead to the transition rate per volume. 

In this paper we concentrate on the calculation of A, we want to find the 
bounce solution numerically. This means to find classical fields Af(t, r) , <E>(i, r) 
which represent a stationary point of the Euclidean action 

S E = f° dtm^(T kili + V&) . (2.9) 

At time t = — oo the fields form a metastable ground state with integer iVcs 
and 7km = 0; at t = to the system reaches the turning point at the other 
side of the potential valley (which is — considering the tunneling process 
- the escape point at the other side of the barrier). Here again = 0, 
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and iVcs takes a non-integer value belonging to the next topological sector. 
The bounce lasts infinitely long because it starts from a ground state with 
5Vp Qt /5Af = 6Vp 0t /8$ = 0. Consequently, to is arbitrary, this corresponds to 
the translational invariance in time mentioned above. 

As discussed previously, the choice of the temporal gauge still leaves the 
freedom of time independent transformations. Now we fix the gauge completely 
by demanding that at t = — oo the fields start from the trivial vacuum Af = 
0, * = (°) with iV C s = (and hence V£ ot = T kin = 0). Besides, we choose 
[i < so that the bounce moves from Nqs = to the topological sector with 
N cs = 1. 

So far, we deal with 13 real functions depending on time t and space r, nine 
from the gauge field and four from the complex Higgs doublet But we 
expect that our bounce solution possesses higher symmetries than completely 
arbitrary fields. Therefore we restrict our ansatz to fields having the spherical 
symmetry of the sphaleron and the AKY configurations, which describe the 
minimal energy barrier: 



$(i,r) 



1 - A(t,r) 



+ (S ai - n a rii) 



B(t,r) C(t,r) 
h n a n,i 



H(t,r) +iG(t,r)n- 



(2.10) 



with r = |r| , n = r/r . This reduces the effort to five real functions depending 
on t and r. The free choice of the origin corresponds to the spatial translational 
invariance discussed after eq. ( |2.8[) . Gauge transformations within this ansatz 
are given by 



U(r) = exp n • 0iP{r) 



they transform the fields as 



(2.11) 



A(t,r 
B(t,r 
C(t,r 
H(t,r 
G(t,r 



A(t, r) cos 2P(r) - B(t, r) sin 2P(r) 
B(t, r) cos 2P(r) + A(t, r) sin 2P(r) 
C(t,r) + 2rP'(r), 
H(t, r) cos P{r) — G(t, r) sin P{r) , 
G(t, r) cos P{r) + H (t, r) sin P(r) . 



(2.12) 



We will use this transformation later to adjust the numerical solutions to our 
choice of the gauge. 
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As mentioned above, we choose a gauge which yields the trivial vacuum at 
t = — oo. In this gauge the fields Af, $ of the bounce solution are continu- 
ous and differentiable everywhere and have finite potential and kinetic energy, 
because this is true for our starting point and will not be changed during the 
evolution governed by the Euclidean equations of motion. This requires the 
following behavior of the radial functions at r = : 

A(t,r) = l + a 2 {t)r 2 + 0{r 3 ), 

B(t,r) = &x(t)r + 0(r 3 ), 

C(t,r) = 6i(t)r + C(r 3 ), (2.13) 

H(t,r) = h (t) + h 2 (t)r 2 + O(r 3 ), 

G(t,r) = gi(t)r + g 2 (t)r 2 + 0(r 3 ). 



The numerical determination of the bounce is performed by finding a sta- 
tionary point of the Euclidean action directly, without using the equations of 



motion. In our spherical ansatz (2.10), Se is a functional of the five functions 
A, B,C, H,G of eq. ( 2.1 0| ) , depending on radial distance r and time t. The 
bounce has infinite extension as well in space as in time, but we can introduce 
new variables x and u which cover only finite intervals, for example 



r{x) = X r arctan (? x) and t(u) = Xt arctan (| u) , 



(2.14) 



and new profile functions depending on x and u. Using ansatz ( |2.10| ) and the 
substitution ( [2.14 ) we get 

Sinst 



Se 



2tt 



du 



dx — 

UJlf 

2 



'(A 2 + B 2 + ^ + 2r 2 {H 2 + G 2 



Y 



+ (yA' + Sf) + (W - 4£) + 2r 2 (^H' + % 
+ 2r 2 UG -^)\ + b2 - !) 2 + H 2 ((A-1) 2 + B 



2r J 
2 +B 2 



4BGH + Wr 2 (H 2 + G 2 - I) 2 



+ c 2 ((A + r 

+ 2p(^(A 2 + B 2 -1)+ ipBA' - ip(A - l)B' 



with 



/'c 



ip = <p(x) 



dry 1 

dx J 



uj(u) 



cky 1 

du J 



(2.15) 



(2.16) 



and the dot and prime mean and 4-, respectively. In the next section we 
describe how we find stationary points of the functional ( p. 15 ) numerically. 



10 



3 Numerical determination of the bounce 
trajectory 

Our way to find a stationary point of Se is the use of a relaxation method which 



was discussed by Adler and Piran in great detail pq |. The functional (2.15) is 
put on a two-dimensional grid of size (n u + 1) x [n x + 1). We distinguish the 
full node grid with points 

+ i Au, x min + j Ax) , . (3.1) 

\ J — u, . . . , n x j 

and the half node grid with points 

x. i) = (u min + (i + ±)Au, x min + (j + ±)Ax) , ( . '•">"•" 
2 J+2 \J = 0, . . . ,n x - 1 

(3.2 



with 



A^-max ^min , A ^max \ 
u = and Ax = . (-J'-jJ 

Tin Tlx 



For the substitution (|2.14[ ) the values are e.g. u m \ n = — 1 , u max = (for t = 0), 
£min = , x max = 1 . The five profile functions are put on the full node grid, 
the radial distance r and the derivatives u and ip on the half node grid. We 
introduce the notation 

A] = A(ui, Xj ) , r j+ i = r(x j+ i) (3.4) 

and accordingly for the other functions. The use of the half node grid prevents 
the calculation of expressions like r or 1/r at the boundaries r = oo,0. 



There are many different ways to put an integral as ( 2.15 ) on a grid, so that 
in the limit n u ,n x — > oo the original functional is restored again. Therefore 
it is rather important to choose a discretization that is suitable for numerical 
treatment. Basically, we followed the suggestions of |l(| here, but to fix the 
details we had to try and compare different ansatzes. For example, the property 
of some terms to vanish at the origin and cancel the 1/r divergence must not 
be lost by the discretization. It showed that terms as (ipA' + ^-) should be 
discretized before being squared, and that instead of C(t,r) the function 

DM = 9M (3.5) 
r(x) 

should be used, i.e. put on the full node grid as D 1 -. We use the following 
discretization for Se, eq. ( 2,15| ), written here as sum over contributions from 
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grid cells centered at half node points (it 4 



grid 



iw^ 1 "^ 1 A« A.i 
2vr 



E E 



i=0 j=0 



2Au 2 



^}) 2 



+ - Aj +1 ) a + - + - ^ + i) 



>i+l 



G}) 2 



+ - G* +1 ) 2 + ipf 1 - D*) 2 + \(Df + \ - D} +1 ) 2 )) 



+ 



Ax 



V Ax v ~' +1 



(Bj +1 -Bi)-^(A5D} + A5 +1 D} +1 



+ 2r U (^ {H ^ Gj+i^+i)) 5 



+ 



3+jV Ax 
1 



(^ + <-l) 2 + (4 +1 2 + i?j +1 2 -l)- 



+ 1 (flf ((AJ - I) 2 + #f ) + Fj +1 2 - I) 2 + B] + *) 

+ Gf ((A) + l) 2 + Ef) + G) +1 \{A) +1 + I) 2 + i?j +1 2 ) 



/I R 1 /I R l z^ 14 IT* 



4 i^+i ( + - 1) 2 + ;//;., 2 + ^.r - 1 r 

+ p( D*(A? + Bf-1) + D) +1 {A) +l 2 + /i)., 1 



1) 



(3.6) 



+ 2 



i -> « + 1 



Typically, we used grids of the size n M = n x = 80, in which case the profile 
functions are represented by 5 • 81 • 81 = 32805 points. To find a stationary 
point of Sg ld now means to find a configuration which fulfills the equations 

n E . = , E . = , E . = , n E . = , n E . = 3.7) 
dA) dB) 3D) dH) dG) v ' 

for all i and j. 
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The solution is found iteratively until the configuration satisfies eq. 
sufficiently well. In each step, only one single number Aj , Bj , Dj , iJj or 
G l j for certain i,j is changed, all others are held constant. "Sweeping" over 
the grid, we modify the field parameters one after another, where the order is 
only of minor importance. (We changed the five functions for given i,j and 
varied j for fixed i.) The change of e.g. is governed only by the according 
equation dSg ld /dAj = 0, all other equations are not taken into account. We 
perform the first step of a Newtonian algorithm that would converge against 
the solution of that equation, i.e. we change the field parameter as 

4 _ AS-.**?™,, (,8, 
3 3 d*Sf A /dAf 

where in general for the damping parameter we choose K = 1. The necessary 
partial derivatives in eq. (3.8) must be calculated from ( p.6[ ), we do not write 
them explicitly here. 

Unfortunately, the bounce solution is not a local minimum, but only a saddle 
point of the Euclidean action: The problem is that because of the term //iVcs 
the Euclidean potential —Vp Qt is not bound from above, it can take positive 
values. Due to energy conservation the bounce itself between iVcs = and 1 
cannot have positive potential energy — Vp Qt at any time, its total energy E tot = 
^kin — Vpot i s constant and zero (Fig. 1). But in its vicinity, one can construct 
paths which have positive potential —Vp Qt for some time and which give a 
lower action than the bounce. We found that unrestricted sweeps according 
to Q3,8| ) always lead to configurations of that kind. Once there, the system 
quickly evolves to enormously high winding numbers iVcs and an unlimitedly 
decreasing negative action. Thus, in order to avoid this instability, we have 
to prevent the system from acquiring positive potential —Vp Qt . We did so by 
choosing an initial configuration with non-positive potential, and then rejecting 
all steps ( |3.8| ) which would yield a positive —Vp 0t . If a step is not accepted, it is 
tried again with n divided by 2; this is repeated up to five times before the step is 
completely rejected for the current sweep. Unfortunately, this method requires 
to calculate V£ ot (or actually the change of Vp 0t ) for the current time slice after 
each single step, which slows down the algorithm considerably. There are more 
sophisticated ways to take into account invariances like energy conservation 
[p7sf | , but in our case the simple remedy proved to be the most effective. Within 
the restriction — Vp Qt < the bounce locally minimizes the action. Numerically, 
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we find indeed that the action decreases monotonically and converges to some 
limit. 

As starting configuration we usually take an instanton-like configuration, 

rt sin (3 + r 2 cos (3 



A(t,r 
B(t,r 
D(t,r 
H(t,r 
G(t,r 
with 



cos [3 — '. 

— sin (3 - 
A 2 



+ 



2 v + 1 2 + a 2 
t 



r 2 + f 2 + X 2 ' 

rt cos (3 — r 2 sin (3 
r 2 + t 2 + A 2 
2rt 



(3.9) 



/3(t,r) 



r 2 + t 2 + A 2 
) (l + cos 

7TT 



7rr 



Vr 2 + A 2 



t 2 + A 2 
2r 



sin 



Vr 2 + A 2 ' 



arctan • 



t 



+ 



Vr 2 + A 2 ^ Vr 2 + A 2 

The size parameter A is chosen between 2 and 4. This configuration has po- 
tential —Vp Qt = at t = — oo and — Vp Q% = 
boundary n max of the grid such that at t(u 
line. Hence for -u m i n < u < u max we have 

According to our gauge fixing and to eq. ( |2.13[ ) the fields are fixed at some 
boundaries of our grid. For u = u m in (t = — oo) or x = x max (r = oo) we have 
A = H = 1 , B = D = G = , and for the origin x = x m [ n (r = 0) we know 
A = 1 , B = G = . Accordingly, the following field parameters are held fixed: 



> at t = oo ; we choose the 
max ) the potential crosses the zero 
"<ot<0. 



A] 


= 4 

Tlx 


= Ai 


= 1 


2$ 


= B L 


= bi 


= 


G) 


- U n x 


— ^0 


= 


B] 


= D i 


= 0, 




H] 


= Hi. 


= 1. 





(3.10) 



Unfortunately, if one only sweeps over the grid using ( |3.8| ) , the outcome will 
not be a reasonable bounce solution, but a rather discontinuous and unusable 
configuration. Therefore the computational procedure cannot be run from start 
to end automatically, but requires the controlling and regulating of the user once 
in a while, which makes the work a tedious and long lasting one. The following 
problems arise: 

The Euclidean action ( 2.15| ) contains no term D' = dD/dx . Hence, adjacent 
field parameters Dj and D l j +1 are only weakly coupled and the D field easily 
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ceases to be smooth. Especially close to the origin the field is rather unstable, 



therefore we do not use (3.8) to fix D l and D\ , but adjust them to a linear 
extrapolation through D\ and D\ . Nevertheless, in order to keep the fields 
reasonably smooth, we are forced to smooth them by an averaging procedure 
from time to time, where each field parameter is replaced by the average value of 
itself and some of its neighbors, with weight factors according to their distance. 
Certainly this disturbs the minimization algorithm and generally leads to a 
higher action again. But after some further sweeps the action is down to its 
earlier value again, and the fields are smoother now. 

Another problem is the following: Close to u = it m i n (t = — oo) the factor 
<jj = du/dt which governs the kinetic energy terms is rather small so that 
adjacent time slices are only weakly coupled. Therefore, the fact that we fixed 
the fields to the trivial vacuum at the u = u m in boundary hardly influences the 
configuration at larger times. Instead, we usually see that going from large u 
towards u m \ n the configuration continuously approaches a non-trivial vacuum 
state, and then close to u m \ n the fields show a discontinuous step from this non- 
trivial vacuum to the trivial one. We get rid of this step by performing a gauge 
transformation of the kind ( |2.12 ) which converts the non-trivial vacuum at the 
edge of the step to the trivial one. This gauge transformation is applied to all 
time slices except the ones close to u m i n where the fields are already trivial. 

Apart from the above manipulations which are necessary to keep the con- 
figuration in an acceptable shape, we also have the possibility to accelerate the 
convergence of the algorithm. The solution has to obey the energy conservation 
law £"tot = ^kin — Vpot = 0- Given some arbitrary configuration for which this 
is not the case, one can find a different time parameterization such that energy 
conservation is fulfilled. Practically, we leave the values t(v,i) of the times at 
the grid nodes fixed and determine the fields of the reparametrized configura- 
tion at those times t(ui) by interpolation. By this operation one can gain a 
considerable decrease of the action without performing minimization sweeps. 

Finally we remark that the grid size is not completely fixed, but it is adapted 
to the status of the minimization. Usually we start with a size of 41 x 41 nodes, 
and only when we are already close to the solution we double the grid to 81 x 81 
points. Moreover, we observe that during the sweeps the bounce, especially 
its escape point where the potential reaches zero again, slowly moves towards 
smaller times. In principle, the time scale is arbitrary, as discussed above, but 
for numerical reasons a given configuration has slightly lower action when it is 
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shifted to smaller times. We find that the energy at a few points next to ii max 
becomes zero, so that the escape point does not longer coincide with u max - 
If the number of these points gets too large, we throw away the part of the 
configuration beyond the escape point, which results in a loW6r Umax 

and Hu- 
ll n u becomes too small, we double the grid in the time dimension only, i.e. we 
double n u but leave n x . For the final configuration we adjust the time scale 
such that the escape point is set on the origin, i.e. to = 0. 

Practically, the bounce trajectory is found by switching between the min- 
imization sweeps and one of the manipulations described above. If any and 
which of the regulations should be performed has to be decided by looking at 
the actual field configuration. At first sight this procedure seems to be subjec- 
tive and non-reproducible, but let us remark that before we stop the program, 
the last manipulation is always followed by at least 500 sweeps. Moreover we 
have always checked that the final configuration fulfills the equations of motion 
with excellent accuracy. 



4 The real time evolution after the tunnel- 
ing 

In this section we show how we investigate the behavior of the boson field 
configuration after the barrier penetration and how we analyze the particle 
content of the state. We follow basically the procedure presented in pi]] , but 
instead of the sphaleron we take the escape point of the bounce as starting 
configuration. 

The potential energy of the escape point is larger by than the potential 
energy of the ground state in the corresponding topological sector. Hence, the 
system performs a motion in the real time Minkowski space. The equations of 
motion within the spherical ansatz ( 2.1 0| ) are 



A = ^ A ' + ^-)' -^{A 2 + B 2 -l) + (^ +2p) 

- A{H 2 + G 2 ) + H 2 - G 2 , 
B = (b 1 -^)' -^{A 2 + B 2 -l)-(^ + 2p){A' + ^-) 

-B(H 2 + G 2 ) + 2HG, 
C = %A(B'-ML)- 2 ^(A' + ^-)-C{H 2 + G 2 ) 

+ 2r(HG' - H'G) - ^(A 2 + B 2 - 1) 
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H = l(rH)" + ^(CG + C'Gr + 2CG'r) + ^(AH + BG) 

-^(l + A 2 + B 2 + ^-)-^H(H 2 + G 2 -l), 

G = l(rG)"-^(CH + C'Hr + 2CH'r) + ^{BH-AG) 

-^(l + A 2 + B 2 + q-)-^-G(H 2 + G 2 -l), (4.1) 

where the dot means the derivative with respect to the time t and the prime 
with respect to the radial coordinate r. 

These equations and the condition that at t = the system starts at the 
escape point of the bounce with kinetic energy zero (i.e. ^4(0) = B(0) = C(0) = 
H(0) = G(0) = 0) form a Cauchi problem which is solved by direct integration. 
To do this we discretize the Minkowskian action in the same way as the Eu- 
clidean action, which means we take eq. (|3.6|) and reverse the sign of the terms 
stemming from the potential energy. Moreover we identify t(u) = u instead 
of eq. (|2.14| ), so that lu(u) = 1, but we still keep the relation between r and 
x of ( p.l4p . We obtain discretized equations of motion by deriving the dis- 
cretized Minkowskian action with respect to the field coordinates Ap . . . , G*- . 
These equations are solved for the variables A % ^ , . . . , so that we obtain 
the fields with time index i+1 as a function of those with indices i and 
By iterative application of these equations it is then possible to evaluate the 
propagation of the system from the initial configuration at t = to arbitrary 
positive times. 

The grid of the discretization can be much more dense here than in the case 
of the Euclidean problem. Usually we take 800 steps per time unit and 3000 
grid nodes in the interval from x = to x = 1. We checked that the results are 
stable with respect to a further increase of these parameters. The total time 
how long we follow the propagation of the fields is typically around 20 to 30 
(in units of m^}). 

As in the Euclidean case, special care has to be taken in order to treat the 
fields close to the origin r = adequately. If the numerical solution does not 
exactly fulfill the expansion given by eq. ( |2.13j ), some terms of the r.h.s. of 
eq. Q4.1| ) become singular. Hence, even a slight deviation from this expansion 
increases rather quickly with t and finally results in a strong divergence of the 
fields close to r = 0. Since small numerical errors will always cause this to 
happen we cannot take the discretized equations of motion close to the origin, 
but we rather impose the behavior of ( |2.13| ) by hand. For about the first 50 
of the 3000 points we determine the fields not by the iteration method but by 
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eq. ( 2.13| ) where the coefficients a%(t) etc. are chosen such that the functions 



A, . . . , H are continuous and differentiable at the matching point between the 
numerical solution and the fit ( 2.1 3j ) . 



We find that for most sets of parameters v and p after some time the fields 
perform small oscillations about some vacuum configuration A(r), . . . , G(r) in 
the topological sector of the escape point (which is the one with iVcs = 1 for 
our choice of gauge). An indication for this behavior is that iVcs> Vp Qt , and 
?kin do not change any more with t, the energies take constant values ±^ 
according to the virial theorem. 

The Fourier analysis of the small oscillations is greatly simplified if the 
vacuum about which the fields are fluctuating is the trivial one. For this reason 
we perform a time independent gauge transformation of the type ( 2.12j ) which 



transforms the configuration A(r), . . . , G(r) into the trivial vacuum. In order 
to determine the configuration A(r), . . . , G(r), we start the propagation of the 
fields with the original escape point and average the fields: 

A{r) = -±- fdt'A{t',r), (4.2) 

t — t s Jt s 

and equivalently for the other profiles. t s is some time where the system already 
performs small oscillations. We checked that A(r), . . . , G(r) are independent 
of t if t is large enough and that they in fact represent a vacuum configuration 
with iVcs = 1- 

The gauge transformation which transforms A(r), . . . , G(r) into the trivial 
vacuum changes the Chern-Simons number by ANqs = — 1- We apply it to 
our starting configuration, the escape point, which hereby gets a Chern-Simons 
number between —1 and and the potential energy V^ ot = |//|. Now we start 
the propagation again with the starting configuration in the new gauge. By 
performing the same averaging process again, we can check that the system 
now indeed fluctuates about the trivial vacuum with excellent accuracy, which 
proves that the gauge invariance is correctly reproduced in our numerics. 

For the following we assume that the system has reached the status where 
it can be described by small fluctuations about the trivial vacuum. We denote 
these fluctuations with small letters: 

A(t,r) = l + a(t,r), B(t,r) = b(t,r) , C(t,r) = c(t,r) , 
H(t,r) = l + h(t,r), G(t,r) = f(t,r), (4.3) 
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and expand the energy to second order in the fluctuations: 

r£ = ^fdr(a 2 + b 2 + ^ + 2r 2 (h 2 + f 2 )), 
< 2) = ^fdr{^ + a> 2 + (b>-tf + a 2 + b 2 + ^ + 4f 2 -4bf 
-2rcf + 2r 2 (h' 2 + f' 2 ) + 2v 2 r 2 h 2 + 2p(a'b - b'a + ^) 

*2 = ^ + <ot (2) - (4-4) 

We found that V£ ot ^ and V^, t coincide up to a deviation of less than 1% for 
large t which indicates that the system has perfectly linearized. The potential 
in second order leads to the linear equations of motion 

a = a"-a(l + ^)+2p(b'-f), 
b = (b' -f)' -b + 2f-2pa> , 
c = f(&'-f)-c + 2r/'-4pf , 
rh = (rh)" — v 2 rh , 

rf = (rf)"-±(rcy+ b -- 2 -f. (4.5) 

In the case p ^ 0, for an arbitrary fixed momentum k the solution is 

«k(t,r) = \{p\{t)-${t))rj x {kr), 

b k (t,r) = ($(t)+(3>[(t)+$(t) + ( k (t))rj (kr) 

+ (/#(*) ~ - \f3 k 2 (t) + C k (t)) rj 2 (kr) , 

c k (t,r) = ^( t )+p k 1 (t)+f3 k 2(t)+C k (t))rj (kr) 

- 2 (/? fc (t) - */3*(t) - \${t) + C k (t)) rj 2 (kr) , 
h k {t,r) = J k (t)j (kr), 

fk(t,r) = -§(*#(*) -JC*(*))ii(*r) (4-6) 



3 
2 

with 



$ {£) = f3i{k) sin^i+a;), i = 0, 1, 2, 7 fc (i) = 7 (fc) sin(m+<5), C (t) = Cl t+c 2 , 

u; 2 = k 2 + 1, uj\ 2 = k 2 + 1 ± 2&p, ft 2 = A; 2 + z/ 2 , (4.7) 

and ji(kr) are the spherical Bessel functions. The phase shifts a^, 5 and the 
q are fixed constants depending on the initial conditions. The zero mode in 
eq. ( |4.6| ) is due to the gauge freedom; since we have fixed the gauge its amplitude 
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C (t) is zero. The frequencies Ui, i = 0, 1, 2 are eigenmodes of free gauge bosons 
(in dimensionful units 1 has to be replaced by m 2 ^ in eq. fl4.7| )), ft represents 
the free Higgs particle eigenstates. 

The coefficients can be evaluated by Fourier transformation: 



l\t) 



2A- 2 



4fc 
"3tt 



dr r 2 jo(kr) (H(t, r) 
dr r 2 j\(kr) G(t, r) 



(4.8) 



2k 2 
9?r 



dr r 



dr r 



j (kr) (2B(t,r) + C(t,r)) 
+ j 2 (kr) (2B(t,r)-2C(t,r)) 



jo(kr) (2B(t,r) + C(t,r)) 

- j 2 (kr) {B(t, r) - C(t, r)) ± 3j a (kr) (A(t, r) - 1) 
To find the amplitudes (3i(k), j(k) we Fourier transform the numerical so- 



lution of the equations of motion (4.1) according to eq. (f4.§|). The resulting 
functions /3^(t), 7 fc (i) should oscillate according to eq. ( |4.7[) when the system 
has settled to small fluctuations about the vacuum. This is the case after some 
time i osc (typically ~ 15 m^) with excellent accuracy. Moreover we checked 



that the two different formulas for [3q (t) in eq. ( [4.8D numerically lead to the same 
result. Fitting the functions /3f (t), 7 fc (t) for t > t osc with the corresponding 
sm(ujit + ai), sm(Qt + 5) of eq. ( f4.7[) yields the amplitudes fli(k), 'y(k). 

One obtains for the total energy of the linearized system in momentum 



space (its coordinate space representation is given by eq. (4.4); the numerical 
values coincide up to a deviation of less than 2%): 



E, 



(2) 
tot 



Ew + Eh 
9ir 2 mw 



= mw 
00 dk 



dkew{k)+mw / dken(k) 



(4(k)^(k)+u^0t(k)+ui(k)0^ 







+ 



r 

47T 2 mw 
~1 F ~ 



k 2 



oo Ah 



(4.9) 



Since energy and particle density are related by e(k) = u>(k)n(k) we can extract 
the total number of particles: 

0^2 roo au . 

r°° dk 



4^ 
~7 



k 2 



n(k)Y(k). 



(4.10) 
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5 Results 



In this section we present the numerical results of our calculation. Our model 
contains two free parameters, the Higgs mass v = mn/fn^r and the chemical 
potential p = fJ,/fJ,ait- We performed the calculations for the values v = 0, 1, 
10, and p = —0.2, —0.4, —0.6, and —0.8. Moreover we investigated the case 
v = 1, p = —0.9. The values for v cover a wide range of Higgs masses, but it 
turns out that most results do not depend too much on v. The choice v = 
is certainly not physical since essential features of the model, like spontaneous 
symmetry breaking, disappear. It should therefore be understood as limiting 
case of small masses. In fact we found that the configurations obtained for 
v = and a mass like v = 0.1 are almost identical. Similarly v = 10 is an 
example for a large Higgs mass. 



5.1 Barrier penetration 

In Fig. 2 we show the potential and kinetic energy as well as the Chern-Simons 
number as functions of the time for the bounce trajectories with chemical po- 
tentials p = —0.2 and —0.8 and the Higgs mass v = 1. The arbitrariness of 
the time scale has been removed by setting the time when the system reaches 
the escape point at the other side of the barrier to t = 0. It can be seen that 
the energy conservation law E tot = T^ n — Vp Qt = is excellently fulfilled. For 
p = —0.8 the energies T^ n and V^ ot are much smaller than for p = —0.2. Also 
the value of A^cs at t = is lower in the case p = —0.8. The reason for this 
behavior is that the barrier between the trivial vacuum and the topological 
sector with Ncs — 1 decreases if \p\ is increased so that less action is neces- 
sary to penetrate it and the system can escape at a configuration with lower 
winding number iVcs- The extension of the bounce in time (and also in space, 
see below), however, increases with \p\. For example, in the case p = —0.2 the 
system moves out of the trivial vacuum significantly at about t ~ —2 while this 
happens for p = —0.8 already at about t ~ —4. The Chern-Simons term in 



the functional Se (eq. ( 2.15 )) lowers the action for large sizes, and its influ- 
ence becomes stronger with increasing \p\. Hence for large \p\ configurations 
with large sizes are favored, while for low \p\ the minimum is taken at a field 
configuration with a small size. In the limit p — > one would even obtain a 
configuration with size zero (see below). 

The field configuration at the escape point of the bounce trajectory at t = 
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-4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0 -0.5 0.0 t 




Fig. 2: The potential energy —Vp 0t , the kinetic energy T^ n , and the 
Chern-Simons number Ncs versus time t (in units of m^) for chemical 
potentials p = /x/// CT it = — 0.2 (solid lines) and —0.8 (dashed lines). The 
Higgs mass is v = mu/rnw = 1- 



is interesting by itself, since subsequent calculations like the investigation of 
the real time behavior of the system require only this configuration rather than 
the complete bounce trajectory. In Fig. 3 we have plotted the profile functions 
at t = again for the two cases p = —0.2 and —0.8, and the Higgs mass v = 1. 
What we saw in Fig. 2 for the time t, we find here for the space coordinate 
r: The deviation of the fields from their values in the trivial vacuum is much 
stronger for p = —0.2, but the region where they deviate is less extended. A 
suitable and accurate analytic fit for the configurations at the escape point is 
provided in the appendix. 
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Fig. 3: The profile functions A, B, C (solid lines), and H and G (dashed 
lines) versus the radial distance r (in units of m^}) of the configuration 
at the escape point of the bounce trajectory for p = /u//x cr it = —0.2 and 
p = —0.8 , and the Higgs mass v = mulfnw = 1- 
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Tab. 1: The action 5e/5i ns t of the bounce trajectory for various 
values of the chemical potential p and the Higgs mass v. 
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2 ) results of two different trajectories with the same action 

Tab. 2: The Chern-Simons number Nq§ of the escape point of the 
bounce trajectory for various values of the chemical potential p and 
the Higgs mass v. 

In Table 1 we give the results for the action Se of the tunneling process in units 
of the action S'inst = ^pr of the instanton in pure gauge theory. In Table 2 we 
show the Chern-Simons number Nq£ of the configuration at the escape point. 
For p = —1 the barrier vanishes so that the tunneling process is reduced to a 
single point in the configuration space, namely the trivial vacuum. Hence in 
this case the action and Nq£ are both 0. For p — > the field configuration which 
minimizes the action tends to an instanton with size zero, so that the bounce 
action equals the instanton action and the configuration at the escape point is a 
vacuum with Nqs = 1- In a pure gauge field theory the solution for p = would 
be an instanton of arbitrary size, but here the scale invariance is destroyed by 
the non-zero vacuum expectation value of the Higgs field, so that for a finite size 



the action would be larger than the instanton action [25|. Hence in the limit 
p->0we obtain a trajectory with size zero. For —1 < p < 0, however, due to 
the Chern-Simons term the minimum of the action is taken by a configuration 
with finite size. The accuracy of the data in Tables 1 and 2 can be estimated 
by increasing the density of the lattice (which is usually of size about 81 x 81) 
and the number of sweeps (usually of the order of 10000). We find that the 
numerical error of the results is around 1%. For each set of parameters v, p we 
performed two independent minimizations, starting from two rather different 
field configurations like e.g. instantons with size A = 2 and A = 4 (see eq. (p.9|)). 
For v = 1, both minimizations always ran towards the same bounce trajectory, 
within the given frame of accuracy. For v = and v = 10, however, the two 
minimizations sometimes produced different trajectories, which have the same 
action, but different Nq£ , and also different behavior of the functions Vp 0t (t), 
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Ncs(t)- Using a suitable averaging procedure, one finds that there exists an 
infinite number of different paths, which all have the same action (up to a 
deviation of 1%). Two conclusions are possible: Either the bounce trajectory 
is not unique, i.e the action has a zero mode, or there is a unique, but very 
shallow minimum. 

Unfortunately, the numerical accuracy of our method does not allow to 
distinguish between those two possibilities, but for the following investigations 
this is quite irrelevant, anyway: In the case of a true zero mode the different 
tunneling trajectories will be taken with exactly the same probability while in 
the case of a shallow minimum the tunneling probabilities are almost the same. 
We will see below that the real time evolution of the fields after the tunneling 
leads to deviating bosonic signatures for two different trajectories. In any case, 
if a barrier penetration happens, both results occur (almost) equally likely, so 
that the accuracy of all results is given by the range of values which we obtain 
for different trajectories with the same action. Certainly, the tunneling rate 
is influenced by the volume of a possible transformation group with invariant 
action or the (low) curvature around the minimum, but these only contribute 
to the prefactor B of eq. ( |2.g| ) which is not discussed in this paper. 

Thus, in Tab. 2 we have given both values for if we obtained two 

different configurations with the same action. We find that both the action 
and Nq£ increase if |p| is decreased or the Higgs mass v is increased. The 
reason is that the barrier becomes wider and higher with decreasing \p\ and 
increasing v so that more action is necessary to penetrate through it, and the 
escape point moves further away from the trivial vacuum. We see that both 
quantities are roughly linearly related to p, and, as anticipated, the dependence 
of the results on the Higgs mass is rather weak. 

The suppression factor e~ 2Sinst ~ 1CP 153 of the tunneling rate for p = 
becomes less strong for \p\ > 0, but significant tunneling amplitudes can only 
be obtained for chemical potentials as high as \p\ > 0.9 (for p = —0.9 and 
v = 1 we obtained Se/Si^ = 0.06 and Nq§ = 0.17). One needs, however, a 
matter density which is about 10 6 times larger than the photon density in the 
early universe at the electroweak phase transition, or about 10 18 nuclear matter 
density, in order to correspond to such a large chemical potential. Presently 
it is not known if a matter density of this order has ever existed in the early 
universe. Even if this is not the case, our results still have some physical signif- 
icance because the tunneling rate may be related to the rate of fermion number 
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violation at high particle energies [14|. Hence, taking our results one might be 
able to deduce the probability to observe such a process at a supercollider. 



5.2 Real time evolution 

Next we will describe how the system evolved in real time after the barrier 
penetration. We solved the equations of motion basically as it was done e.g. in 



, 22 1 . In order to check the numerics, we first used the sphaleron as starting 



configuration, and our results agreed with those of [21, E3|. Then we replaced 
the sphaleron by the configuration at the escape point of the barrier penetration 
process. By performing the gauge transformation with AiVcs = — 1 to ensure 
that the fields fluctuate about the trivial vacuum (see Section 4) the potential 
energy of the starting configuration V^ ot is increased from to \fj,\ (in units of 
^p-mw = ^p- this means it is increased from to 2\p\). The Chern-Simons 
number is lowered by one unit and starts between —1 and 0. 

In Fig. 4 we present the behavior of V^ ot , T^[ Q , and Nqs as a function of 
time t > after the tunneling process for p = —0.6 and p = —0.9. The Higgs 
mass is v = 1. We find that for times £ ^ 10 the behavior of the system in the 
two cases is quite similar: The system starts to move, i.e. the kinetic energy 
increases, while the potential energy decreases. Energy conservation is fulfilled 
very accurately during the whole process. The Chern-Simons number increases 
quickly to values around 0, which means the systems comes close to the trivial 
vacuum. 

In the case p = —0.6 the energy dissipates into small fluctuations about the 

rH ( 
'pot 



vacuum. We checked that the energy VT ot ^ in second order of the fluctuations 



agrees up to a deviation of less than 1% with Vp Dt , moreover we see in Fig. 4 
that potential and kinetic energy both become constant at the value of ^ in ac- 
cordance with the virial theorem. The Chern-Simons number takes a constant 
value slightly below zero. Hence the system has settled to small oscillations 
about the trivial vacuum and will stay in this topological sector forever (apart 
from possible tunneling later on). Later we will analyze the particle content of 
this state. 

In the case p = —0.9, however, we observe a completely different behavior 
for times t > 10. Here iVcs suddenly increases from values around to about 
1, later even to 2. Hence, the system does not stay in the topological sector of 
the trivial vacuum but moves classically over the next barrier to the sector with 



26 




~\ — i — i — i — | — i — i — i — i — | — i — i — i — i — | — i — i — i — i — | — i — i — i — i — | — i — i — i — r 

tO 5 10 15 20 25 

J I I I I I I I I I I I I I I I I I I I I I I I I I I I L 



30 




Fig. 4-' The potential energy Vp Qt , the kinetic energy T^ in (in units of 
87r2 ^ w ), and the Chern-Simons number TVcs versus time t (in units of 
m^/ 1 ) for chemical potentials p = p/p C rit = —0.6 (solid lines) and —0.9 
(dashed lines), and the Higgs mass v = mn/mw = 1- 



iVcs = 1- Here it also stays only for a short period before it moves to the next 
sector with Nqs = 2. This behavior is also demonstrated by the plot of the 
potential energy which shows the successive falls of the system like a cascade 
towards configurations with increasing winding number and decreasing energy. 
Once the first tunneling process has happened, the system moves classically 
over all the following barriers so that the whole fermion matter decays rapidly 
and sets free an enormous amount of energy. As mentioned above, the tunneling 
amplitude is not extremely small for p = —0.9 (10~ 9 instead of 10~ 153 ), but 
in order to generate such a large chemical potential a huge fermion density is 
required (f» 10 18 nuclear matter density). 
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It is a property of the periodic plus linear potential that even at small p, 
the energy barriers become lower than the local minimum at Acs = 0, if one 
goes far enough in iV"cs- Therefore, if the systems tunnels directly to that 
far-away sector (which would require multi-instanton-like bounce solutions), 
the avalanche would probably develop, too. Of course, the multi-instanton 
tunneling probability is even smaller than for a single bounce, but it should grow 
faster with p. It would be interesting to estimate the total decay probability 
as a function of p, with tunneling to different topological sectors summed up. 

Energetically, this behavior of a classical rapid decay is allowed if the top of 
the barrier between the sectors with Nqs = and 1 is lower than the chemical 
potential \p\. For v = 1, this is the case already for |p| <; 0.2, but we found 
that only for \p\ Z 0.9 it actually happens. For chemical potentials between 0.2 
and 0.9 the system could in principle cross the next barrier, but the energy is 
dissipated among the modes of small oscillations and not concentrated on the 
direction to the next minimum so that the system does not find the collective 
path over the barrier. For v = 0, we have found that the avalanche starts 
developing already at \p\ > 0.8. We think we have observed an interesting 
phenomenon of how an exponentially suppressed spontaneous decay triggers 
off a catastrophic avalanche which never stops until the fermion sea, originally 
filled up to the Fermi surface p, is completely "splashed" . 

In Fig. 5 we show for v = 1 and p = —0.6 how the density of the total energy, 
defined by E to t = Ikin + Vp Qt = my/ J °° dr e to t(^)j evolves in time. Our plot is 
similar to the one given in [21, |2^], where one starts with a slightly disturbed 
sphaleron instead of the configuration at the escape point of the bounce. As 
was found in |2l], |22|], the outgoing wave moves with almost the speed of light 
and shows some dispersion, but in our case the dispersion is less strong. For 
the sphaleron, after t = 25 the height of the pulse has decreased to about 
30% of its original value at t = 0, while for the bounce configuration it drops 
only to about 65% of the value at t = 0. 

Next we wish to analyze the particle content of the state after the tunneling. 
This is only possible if the system stays in the topological sector of the trivial 
vacuum and does not move classically over the next barrier. In this case after 
some time (typically ~ 10 m^/) the system has settled to small oscillations 
about the trivial vacuum so that one can perform the Fourier decomposition 
( ^4.8| ). Fig. 6 shows how the total energy (= \p\) is distributed among the Higgs 
and gauge boson modes. Integration of the curves yields the total energy of 
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Fig. 5: The density etot of the total energy versus radial distance r 
for various times t (r and t in units of m^}). The parameters are 
p = jU//Xcrit = —0.6 and v = mn/mw = 1. 



the Higgs (Eh) and gauge (E\y) bosons (see eq. (4.9)); we have found that the 
sum Eh + E\y is equal to \fj,\ up to a deviation of usually less than 2%, which 
is another check of our numerics. 
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1 ^ results of two different trajectories with the same action 
2 ) system does not oscillate about trivial vacuum 



but moves classically to next sector 

Tab. 3: The ratio of the energy of the Higgs bosons to the total 
energy Eh/(E\y + Eh) in percent after the system settles to small 
oscillations about the trivial vacuum. The results are given for var- 
ious values of the chemical potential p and the Higgs mass v. 
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In Tab. 3 we show how much of the energy is taken by the Higgs bosons (in 
percent). This number is generally in the range up to 10%; it increases slightly 
with \p\. If we increase the Higgs mass from v = to u = 1 (and keep p fixed), 
the Higgs particles gain some energy on the expense of the gauge bosons, but 
for v = 10 the share of the Higgs bosons is almost zero. In this last case the 
Higgs bosons are too heavy to be produced at all, for small masses their total 
energy is basically correlated to the individual energy of each particle, i.e. it 
rises with the mass. 

Fig. 6 demonstrates that the spectrum is shifted to larger k when v is 
increased. This effect is particularly strong in the case of the Higgs bosons. For 
v = 10 the energy density takes its maximum at about k « 5 while for v < 1 
it is only at k « 1.5, for v = very light Higgs bosons with momenta around 
k = are produced in a large number. 




Fig. 6: The energy densities e\y and en of the gauge and Higgs bosons 
versus momentum k (in units of my/) for Higgs masses v = mn/mw = 
(solid lines), v = 1 (dashed lines), and v = 10 (dotted lines). The 
chemical potential is p = —0.6. 
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Tab. 4 shows the total particle numbers Nyp and Njj of the gauge and 
Higgs bosons, respectively. In the case v = the determination of Nh is not 
possible because the number density njf(k) is strongly peaked close to k = 
so that the numerical error of the integration is uncontrollable. This is of 
course an artefact of the unphysical choice v = 0, however, also for finite, but 
small masses one would have to take many values in the fc-lattice around 
and perform the integration carefully to get a reasonable result. We find that 
the particle numbers rise with \p\, because the released energy increases and 
allows the production of more particles. The number of Higgs particles is much 
smaller than the number of gauge bosons, again we see that in case of a large 
Higgs mass no Higgs bosons are produced. 
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Tab. 4: The number of gauge bosons Ny/ (upper numbers) and Higgs 
bosons Nh (lower numbers) after the system settles to small oscil- 
lations about the trivial vacuum. The results are given for various 
values of the chemical potential p and the Higgs mass v. 

Finally we comment on the numerical uncertainty of the data in Tabs. 3 and 
4. The error can be estimated by increasing the numerical parameters of the 
Runge-Kutta time integration and the number of k values in the Fourier trans- 
formation. Moreover one can choose different times i osc (see explanation after 
eq. ([4.8])) where we start to fit the amplitudes of the modes in momentum 
space. We find that the error of the data is in general less than 2%. We have 



31 



to keep in mind, however, that the bounce trajectory, and hence the starting 
configuration, is not unique, but in some cases there are different solutions with 
almost the same tunneling probability. These different starting configurations 
yield results for the particle content which can deviate up to 20%, as can be 
seen from the data in Tabs. 3 and 4. Therefore, the probability density for 
the particle content of the final state is spread over a range of numbers about 
±10% around the value given in the tables. 

6 Summary 

In this work we have presented a method to find the bounce trajectory in the 
electroweak theory and calculated the probability for the decay of high density 
fermionic matter. 

The bounce trajectory is obtained by minimization of the Euclidean action 
as a function of the discretized Higgs and gauge boson fields. At each step 
of the procedure the action is regarded as a function of only one parameter, 
i.e. it is minimized with respect to the value of one field profile function at a 
certain point in the lattice while the values of the other fields and at the other 
points are kept fixed. After finishing one step of the minimization one moves 
to the next field or next point until each field at each point of the lattice has 
been considered. Many of such "sweeps" through the lattice (of the order of 
10000) have to be performed until a stable configuration is reached which does 
not change any more if it undergoes further sweeps. From time to time the 
user has to interfere into the process of minimization. The program contains 
several options to manipulate the field configuration, partially in order to keep 
the fields in a continuous and smooth shape, partially in order to accelerate the 
convergence. It has been checked that the final configuration always fulfills the 
Euclidean equations of motion with sufficient accuracy. 

The determination of the bounce has been carried out for several choices 
of the Higgs mass u = m^/m^ and the chemical potential p = p/p CT it of the 
fermionic matter. We find that the action Se of the bounce drops from the 
instanton action Si ns t at p = to zero at \p\ = 1 roughly linearly and depends 
only weakly on v. A similar behavior is found for the Chern-Simons number 
of the escape point of the bounce, iVgg, which decreases from 1 to 0. The 
action Se is the exponent of the tunneling rate which itself is correlated to the 
probability of the fermion number violation at high particle energies. It might 
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therefore be possible to use our results for Se in order to predict the cross 
section of the high energy process. 

For several sets of parameters we found that the bounce solution is not 
unique; instead there exist several solutions with different escape points, but 
with the same action (in the given frame of accuracy). Since after the barrier 
penetration each of these escape points will be taken with the same probability, 
for some results of this work we can only give a range of values instead of a 
definite number. 

After the tunneling process, the bosonic fields can evolve in real time Min- 
kowskian space since they have obtained the energy of the annihilated fermions 
as potential energy. The equations of motion can here be solved by some time 
integration method rather than a minimization of the action. 

For chemical potentials \p\ < 0.8 the system stays in the topological sector 
where it came to after the tunneling and settles to small oscillations about 
the minimum. The oscillations correspond to the radiation of the Higgs and 
gauge bosons, and we have analyzed the particle content of this state by Fourier 
transformation. We find that usually less than 10% of the energy is absorbed 
by the creation of the Higgs bosons, and correspondingly the total number 
of produced gauge bosons is also about 10 times greater than the number of 
Higgs bosons. The results depend strongly on p and partially also on the Higgs 
mass v. For large chemical potentials \p\ > 0.9 the system has enough energy 
and coherence after the tunneling to move classically over the next barriers. 
This corresponds to an avalanche decay of the fermionic matter and to the 
production of an enormous amount of Higgs and gauge bosons. 
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Appendix 



To investigate the real time behavior after tunneling one only needs the field 
configuration at the escape point of the bounce (which we fix at t = 0) rather 
than the functions A(t,r), . . . , G(t,r) in the whole two-dimensional space t,r. 
Hence it is useful to have a parameterization of the functions ^4(0, r), . . . , G(0, r) 
so that one can take them as input for further calculations without having 
the necessity to recalculate the complete bounce trajectory. In this appendix 
we give an analytic fit which matches the numerically determined functions 
A(0,r), G(0,r) very accurately. The potential energy and the Chern- 
Simons number of the fit agree with the corresponding values of the numerical 
configuration up to a deviation of about 1%. Following the real time behavior 
of the fit and the numerical fields we find that the particle numbers Nyy, Nh 
and the energies E\y, Eh are reproduced up to a deviation less than 3%. 

The parameterization is chosen so that it includes the possibility to describe 
both the trivial vacuum with iVcs = and the non-trivial vacuum with Nqs 
1 . For this reason the fit is performed in a gauge where the field D(r) = 
everywhere. We denote the other fields in this gauge by Ao(r), Bo(r), Hq(v), 
Go(r). (Here and in the following the argument t = is dropped.) The gauge 
transformation which transforms ^4o( r )> • • • > Go(r) to the field configuration 
A(r), . . . , G{r) at the escape point of the bounce is described by some function 
P(r) according to eq. (|2.12|) . 

For the functions P(r), D(r), Aq(t), Bo(r), Ho(r), Gq{t) we use the follow- 
ing ansatz: 
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(A.l) 



The procedure how we obtain the parameters in eq. (|A.1|) is the following: To 
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ensure the correct behavior of the fitting functions at r = (eq. ( |2.13[ )) we first 
take do = B'(0) (evaluated by a quadratic fit of B at r = 0), use a suitable 
fitting algorithm for D(r) to determine A^, d 2 , and d 3 , and set 

a = cos(2X d (d + d 2 + 3d 3 )) - 1 , 

6 = sin(2A d (do + d 2 + 3c k)) , 

h = H(0)cos(\ d (d + d 2 + 3d 3 f) - 1, 

go = H(0) sm(\ d (d + d 2 + 3d 3 )) , 

hi = (H(0)f-G'(0))sm[\ d (do + d 2 + 3d 3 )), 

91 = -(H(0)f-G'(0))cos(\ d (do + d 2 + 3d 3 )), (A.2) 

where we have used P(0) = — A^(do + ^2 + 3^). The parameters A a , A&, A/ t , 
A g , are held fixed: 

A a = A b = 0.8, A/, = Aj, = 0.6. (A.3) 

Then we perform a gauge transformation on the fields A(r), . . . , G(r), using 
the function P(r) of eq. ( |A.1| ) to obtain Ao(r), ■ • • , Go(r). These functions are 
fitted to yield the remaining parameters a 2 , (Z3, b 2l 63, h 2 , /13, g 2 , 53, and #4. 
Altogether, our fits contain 12 parameters determined by the fitting algorithm 
plus 7 parameters depending on the three values B'(0), H(0), and G'(0). So in 
total the number of free parameters is 15. 

In Table 5 the results of the parameters are given for several values of the 
chemical potential p and the fixed Higgs mass v = mH/mw = 1- For this 
Higgs mass we always obtain a unique field configuration at the escape point, 
i.e. it does not depend on the initial configuration before the minimization of 
the action starts. 
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Tab. 5: The parameters of the fitting functions for the configuration 



at the escape point of the bounce trajectory according to eq. ( A.l ) 
for several values of the chemical potential p and the Higgs mass 
v = run/raw = 1- 
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